
rm(list=ls(all=TRUE))
library(ggplot2)
library(gridExtra)
library(car)
library(plm)
library(gam)
library(plyr)
library(doBy)
library(foreign)
library(texreg)
library(xtable)
library(MASS)
library(GGally)
library(zoo)
library(lmtest)
library(plm)
library(lattice)
library(ebal)
library(DataCombine)
library(rgenoud)
library(sandwich)
library(Hmisc)
library(dplyr)
library(tidyr)
library(lubridate)
library(zoo)
require(texreg)
require(lmtest)
require(multiwayvcov)
library(VIM)
library(rpart)

setwd("/Users/irismalone/Dropbox/AGDPapers/IS-ISSS Paper/Data/")

##
#moral hazard -
dfsumm=read.csv("/Users/irismalone/Dropbox/Dissertation/Data/TORG/January/fullsamplebycountryv2.csv")
names(dfsumm)
summary(dfsumm$noact)

dfsumm=read.csv("/Users/irismalone/Dropbox/Dissertation/Data/TORG/January/fullsamplecountryyearv14.csv")
names(dfsumm)
summary(dfsumm$noactviol)

##
df = read.csv("/Users/irismalone/Dropbox/AGDPapers/IS-ISSS Paper/Data/dyadicextsupp_fullv2.csv")
df = subset(df, df$torgid != 28)
df$anystatesup=as.factor(df$statesup)
df$year = ifelse(is.na(df$yrfirstattack)==TRUE, df$yrform, df$yrfirstattack)
df$coldwar = ifelse(df$year < 1991, 1, 0)
df$coldwar[is.na(df$year) == TRUE & df$yrform < 1991] = 1
df$post911 = ifelse(df$year > 2001, 1, 0)
df$post911[is.na(df$year) == TRUE & df$yrform > 2001] = 1
df$usasponsor = ifelse(df$ccode2==2, 1, 0)
df$russponsor = ifelse(df$ccode2==365, 1, 0)

df$aim2 = df$aim2-1
df=with(df, data.frame(ccode1, ccode2, torgid,  anystatesup, statesup, coldwar, post911, usasponsor, russponsor, 
                       majpow1, majpow2, pol_rel, contig, distance))
dim(df)


dftorg = read.csv("/users/irismalone/Dropbox/Dissertation/Data/TORG/January/dfsamplecleanfullv4.csv")
dftorg$year = ifelse(is.na(dftorg$yrfirstattack)==TRUE, dftorg$yrform, dftorg$yrfirstattack)
table(dftorg$statesup[dftorg$polopp==1])
dftorg = subset(dftorg, dftorg$polopp==1)
table(length(unique(dftorg$torgid[dftorg$statesup==1])))
table(length(unique(dftorg$ccode[dftorg$statesup==1])))

table(length(unique(dftorg$ccode)))
table(length(unique(dftorg$torgid)))
dftorg$numstatesponsors[dftorg$numstatesponsors==-9] = NA
dftorg = dftorg %>% select(-numnonstateties.1, -numstatesponsors.1, -numbases.1)


dftorg$ccode1 = dftorg$ccode
df = merge(df, dftorg, by=c("ccode1", "torgid"), all=T)
head(df)

df$left[df$torgid==24] = 0
df$leftideational = ifelse(df$left==1 & (df$ccode2== 40 | (df$ccode2==365 & df$year <1992) | df$ccode2 == 710 | (df$ccode2==530 & df$year <1987) | (df$ccode2==514 & df$year <1990) ), 1, 0)
df$shiaideational = ifelse(df$religious==1 & df$torg != "T28" & (df$ccode2==630 | df$ccode2==652), 1, 0)
df$sunniideational = ifelse(df$religious==1 & df$torg != "T208" & (df$ccode2==651 | df$ccode2 ==645 | df$ccode2==690 | df$ccode2 == 694 | df$ccode2 == 670 | df$ccode2==770 | df$ccode2==700 | df$ccode2==625), 1, 0)
df$sharedrelig = df$shiaideational+df$sunniideational



df2 = read.csv("/Users/irismalone/Dropbox/AGDPapers/IS-ISSS Paper/Data/eth_shared_ties2.csv")
table(df2$sharedethnicity)
df2 = with(df2, data.frame(ccode1, ccode2, torgid, sharedethnicity, ccode2epr))

df = merge(df, df2, by=c("ccode1", "ccode2", "torgid"), all=T)
df = subset(df, is.na(df$ccode2) == FALSE)
df = subset(df, is.na(df$ccode1) == FALSE)
df$ccode2epr[df$ccode2epr=="<NA>"] = NA
df$sharedethnicity[is.na(df$ccode2epr)==TRUE] = 0
df$sharedethnicity[df$ethmem==0 & is.na(df$sharedethnicity)==TRUE] = 0
df$sharedethnicity[df$ethmem==0] = 0

dim(df)
table(is.na(df$sharedethnicity), df$ethmem)

require(foreign)
#colgan = read.csv("/Users/irismalone/Dropbox/AGDPapers/IS-ISSS Paper/Data/colgan_ideo.csv")
#colgan$ccode2 = colgan$ccode
#colgan$yrfirstattack = colgan$year
#df = merge(df, colgan, by=c("ccode2", "year"), all.x=T)

table(df$pol_rel, df$statesup.x)
table(df$anystatesup)

dfislam = read.csv("/users/irismalone/Downloads/Islamic.csv")
dfislam = with(dfislam, data.frame(torgid, shia, sunni, salafistjihadist))
dfislam$ccode1 = dfislam$ccode
df = merge(df, dfislam, by=c("torgid"), all.x=T)
df$shia[df$islam==0] = 0
df$sunni[df$islam==0] = 0
df$salafistjihadist[df$islam==0] = 0

df$sunni[df$islam==1 & df$ideo2 == "Islamic-NA" & df$shia==0] = 1
df$shia[df$islam==1 & df$sunni==0] = 1
df$sunni[df$islam==1 & is.na(df$shia)==TRUE] = 1
df$shia[df$islam==1 & df$sunni==1] = 0
df$sunni[df$islam==1 & df$shia==1] = 0

#df$sunni[df$islam==1 & df$ideo2 == "Islamic-NA" & df$shia==0] = 1

table(df$sunni, df$shia)


df$leftideational = ifelse(df$left==1 & (df$ccode2== 40 | (df$ccode2==365 & df$year <1992)  | df$ccode2 == 731 | df$ccode2 == 710 | (df$ccode2==530& df$year <1987) | (df$ccode2==514 & df$year <1990) ), 1, 0)
df$leftideational[df$left==0] = 0
table(df$left, df$leftideational)
table(is.na(df$leftideational))

df$shiaideational = ifelse(df$shia==1 & df$torg != "T28" & (df$ccode2==630 | df$ccode2==652), 1, 0)
df$shiaideational[df$shia==0] = 0
df$sunniideational = ifelse(df$sunni==1 & df$torg != "T208" & (df$ccode2==651 | df$ccode2==600| df$ccode2==663 | df$ccode2==620 | 
                                                                 df$ccode2==616 | df$ccode2==620 |df$ccode2 ==645 | df$ccode2==690 | 
                                                                 df$ccode2 == 694 | df$ccode2 == 670 | df$ccode2==770 | 
                                                                 df$ccode2==700 | df$ccode2==625 | df$ccode2==696 | df$ccode2==692 | df$ccode2==820 | df$ccode2==850), 1, 0)
df$sunniideational[df$sunni==0] = 0
df$sharedrelig = df$shiaideational+df$sunniideational
df$ccode2epr[df$ccode2epr=="jr partner" | df$ccode2epr=="sr partner" | df$ccode2epr=="senior partner"] = "some power"
df$sharedethnicpower = ifelse(df$sharedethnicity == 1 & (df$ccode2epr == "somepower" | df$ccode2epr == "dominant" | df$ccode2epr == "monopoly" ), 1, 0)
df$sharedsomepower = ifelse(df$sharedethnicity == 1 & df$ccode2epr == "somepower", 1, 0)
df$shareddominant = ifelse(df$sharedethnicity == 1 & (df$ccode2epr == "dominant" | df$ccode2epr=="monopoly"), 1, 0)
df$sharednoethnicpower = ifelse(df$sharedethnicity == 1 & (df$ccode2epr == "autonomous" | df$ccode2epr == "discriminated" | df$ccode2epr == "irrelevant" | df$ccode2epr == "powerless" ), 1, 0)
table(df$sharedrelig)
table(df$shiaideational)
table(df$sunniideational)
table(df$leftideational)

df$sharedcoideo = df$leftideational+df$sharedrelig+df$sharedethnicpower
table(is.na(df$leftideational))
df$sharedcoideo[df$sharedcoideo >0] = 1

df3=with(df, data.frame(ccode1, ccode2, year, torgid, group, ccode2epr, polopp, anystatesup, left, islam, shia, sunni, leftideational, sharedethnicity, sharedethnicpower, sharedrelig, sharedcoideo ))
#df = merge(df, dfnoarmed, by=c("ccode1", "startyear"), all=T)

#######
#summary statistics


df$lndist = log(df$distance)
library(mice)

table(df$sharedrelig)
table(df$sharedethnicity)

table(df$sharedethnicity, df$statesup.x)
table(df$sharedrelig, df$statesup.x)
df$coldwar = ifelse(df$year < 1991, 1, 0)
df$post911 = ifelse(df$year > 2001, 1, 0)

require(stargazer)
df$sharedcoideo = df$leftideational+df$sharedrelig+df$sharedethnicpower
#df$sharedcoideo = df$leftideational+df$sharedrelig
df$sharedcoideo[df$sharedcoideo >0] = 1

df2 = subset(df, df$pol_rel==1 & df$torgid!=28 & df$polopp==1)

df2$anystatesup=as.numeric(df2$anystatesup)-1

dfsummstats = with(df2, data.frame(anystatesup, sharedethnicpower, sharedcoideo, leftideational, sharedrelig, shiaideational, sunniideational))
stargazer(dfsummstats, median = T, omit.summary.stat = c("p25", "p75"))
table(df2$sharedrelig, df2$statesup.x)
t.test(df2$statesup.x ~ df2$sharedcoideo)

t.test(df2$statesup.x ~ df2$sharedethnicity)
t.test(df2$statesup.x ~ df2$sharedethnicpower)
t.test(df2$statesup.x ~ df2$sharedrelig)

table(df$sharedethnicity, df$statesup.x)

table(df2$sharedcoideo[df2$coldwar==1], df2$statesup.x[df2$coldwar==1])
t.test(df2$statesup.x[df2$coldwar==1] ~ df2$sharedcoideo[df2$coldwar==1])
t.test(df2$statesup.x[df2$coldwar==0] ~ df2$leftideational[df2$coldwar==0])
require(dplyr)
library(tictoc)
library(peacesciencer)
tic()
create_dyadyears(directed = TRUE, mry = FALSE) %>%
  #filter_prd() %>%
  add_gml_mids(keep = NULL) %>%
  add_peace_years() %>%
  add_nmc() %>%
  add_democracy() %>%
  add_cow_majors() %>%
  add_cow_alliance() %>%
  add_contiguity() %>%
  add_strategic_rivalries() %>%
  add_sdp_gdp() -> Data

require(dplyr)
Data %>%
  mutate(rivalry = ifelse(is.na(rivalryno)==FALSE, 1, 0)) %>%
  mutate(landcontig = ifelse(conttype == 1, 1, 0)) %>%
  mutate(cowmajdyad = ifelse(cowmaj1 == 1 | cowmaj2 == 1, 1, 0)) %>%
  # Create estimate of militarization as milper/tpop
  # Then make a weak-link
  mutate(milit1 = milper1/tpop1,
         milit2 = milper2/tpop2,
         minmilit = ifelse(milit1 > milit2,
                           milit2, milit1)) %>%
  # create CINC proportion (sponsor over target. if negative then sponsor weaker than target)
  mutate(cincprop = 
           cinc2/cinc1) %>%
  # create weak-link specification using Quick UDS data
  mutate(mindemest = ifelse(xm_qudsest1 > xm_qudsest2,
                            xm_qudsest2, xm_qudsest1)) %>%
  # Create "weak-link" measure of jointly advanced economies
  mutate(minwbgdppc = ifelse(wbgdppc2011est1 > wbgdppc2011est2,
                             wbgdppc2011est2, wbgdppc2011est1)) -> Data

dim(Data)
head(Data)
table(Data$rivalry)

#Data$dyadid = Data$ccode2*1000 +Data$ccode1
Data = subset(Data, Data$year >= 1945)

head(Data)
head(df)
df = merge(df, Data, by=c("ccode1", "ccode2", "year"), all=T)


require(lmtest)
require(multiwayvcov)

dfsumm=read.csv("/Users/irismalone/Dropbox/Dissertation/Data/TORG/January/fullsamplecountryyearv13.csv")

dfsumm = dfsumm %>%                            
  group_by(ccode) %>%
  dplyr::mutate(l1noactviol = lag(noactviol, n = 1, default = NA)) %>%
  dplyr::mutate(l1noviolsup = lag(noviolsup, n = 1, default = NA)) %>%
  dplyr::mutate(l1noleftsup = lag(noleftsup, n = 1, default = NA)) %>%
  dplyr::mutate(l1noethnosup = lag(noethnosup, n = 1, default = NA)) %>%
  dplyr::mutate(l1noislamsup = lag(noislamsup, n = 1, default = NA)) %>%
  dplyr::mutate(l1noideosup = lag(noideosup, n = 1, default = NA)) %>%
  dplyr::mutate(l1noactislam = lag(noactviolislam, n = 1, default = NA)) %>%
  dplyr::mutate(l1noactleft = lag(noactviolleft, n = 1, default = NA)) %>%
  dplyr::mutate(l1noactethno = lag(noactviolethno, n = 1, default = NA)) %>%
  dplyr::mutate(l1noactiveinsurg = lag(noactiveinsurg, n = 1, default = NA)) 



dfsumm$ccode1 = dfsumm$ccode
dfsumm2 = with(dfsumm, data.frame(ccode1, year,  l1noactviol, l1noviolsup, l1noleftsup,
                                  l1noislamsup, l1noethnosup, l1noideosup, l1noactislam, l1noactleft, l1noactethno, l1noactiveinsurg))
df = merge(df, dfsumm2, by=c("ccode1", "year"), all.x=T)


#df2 = subset(df, df$pol_rel==1 & df$torgid!=28 & df$polopp==1)
df2 = subset(df, df$torgid!=28 & df$polopp==1)

df2$anystatesup=as.numeric(df2$anystatesup)-1
head(df2$dyadid)

df2$numnonstateties[df2$numnonstateties == -9 & df2$nsasup == 1] = 1
df2$numnonstateties[df2$numnonstateties == -9 & df2$nsasup == 0] = NA

df2$sharedcoideo2 = df2$leftideational+df2$sharedrelig
df2$sharedcoideo2[df2$sharedcoideo2 > 1] = 1

write.csv(df2, "replication_ii_fragment3.csv")

write.csv(df2, "externalsupport_robustnessdatav2.csv")


############################
## SUMM STATS
###########################
prop.table(table(df2$sharedcoideo, df2$anystatesup),1)
prop.table(table(df2$sharedcoideo, df2$anystatesup),2)
prop.table(table(df2$sharedcoideo[df2$civwar==1], df2$anystatesup[df2$civwar==1]),2)

prop.table(table(df2$sharedcoideo, df2$anystatesup),2)

prop.table(table(df2$sharedcoideo2, df2$anystatesup),1)
prop.table(table(df2$sharedcoideo2, df2$anystatesup),2)

prop.table(table(df2$l1noideosup, df2$anystatesup),1)

#df2$dyadid = df2$dyadid.x

xtable(table(df2$pol_rel, df2$anystatesup))
prop.table(table(df2$pol_rel, df2$anystatesup), 1)


table(df2$numnonstateties)
table(df2$l1noviolsup)



################
# Hypothesis 1 
################

df2=read.csv("/Users/irismalone/Dropbox/AGDPapers/IS-ISSS Paper/Data/v4/externalsupport_robustnessdata.csv")

df2$l1noviolsupl = log(df2$l1noviolsup+1)
df2$l1nostatesupl = log(df2$numstatesponsors+1)
table(length(unique(df2$torgid)))

table(df2$numstatesponsors)
table(df2$numstatesponsors[df2$anystatesup==1])
length(unique(df2$torgid[df2$anystatesup==1]))
length(unique(df2$torgid[df2$numstatesponsors==1]))
length(unique(df$torgid[df$numstatesponsors==1]))
length(unique(df2$torgid[df2$anystatesup==1]))
unique(df2$torgid[df2$numstatesponsors==8])

cor(df2$numstatesponsors, df2$l1noviolsup, use="complete")

#table 1
reg.m0 = glm((anystatesup) ~ l1noviolsupl + sharedcoideo, data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m0, df2$dyadid)
a = coeftest(reg.m0, vcov_firm)
a

reg.m1 = glm((anystatesup) ~ l1noviolsupl*sharedcoideo, data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m1, df2$dyadid)
b = coeftest(reg.m1, vcov_firm)
b

reg.m2 = glm((anystatesup) ~ l1noviolsupl*sharedcoideo  + rivalry + minwbgdppc + mindemest, data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m2, df2$dyadid)
c = coeftest(reg.m2, vcov_firm)
c

reg.m3 = glm((anystatesup) ~ l1noviolsupl*sharedcoideo + rivalry + minwbgdppc + mindemest + factor(ccode1), data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m3, df2$dyadid)
d = coeftest(reg.m3, vcov_firm)
d


texreg(list(a, b, c, d), omit.coef="ccode1", stars=c(0.1, 0.05, 0.01))
texreg(list(reg.m0, reg.m1, reg.m2, reg.m3))



#############
# predicted probability values

reg.m1 = glm((anystatesup) ~ l1noviolsupl*sharedcoideo + minwbgdppc + mindemest + rivalry, data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m1, df2$dyadid)
b = coeftest(reg.m1, vcov_firm)
b
require(ggplot2)
require(ggeffects)

p_1000 <- ggpredict(reg.m1, c("l1noviolsupl [meansd]", "sharedcoideo"), vcov.fun = "vcovCR", vcov.type = "CR0", 
                    vcov.args = list(cluster = df2$dyadid))
p_1000

p_1000 <- ggpredict(reg.m1, c("l1noviolsupl [0, 0.5, 0.69, 1, 1.5, 2, 2.5, 3, 3.5]", "sharedcoideo"), vcov.fun = "vcovCR", vcov.type = "CR0", 
                    vcov.args = list(cluster = df2$dyadid))
p_1000

a = plot(p_1000, colors="bw") + 
  theme_bw() + ylab("Predicted Pr(External Support)") + xlab(expression("Ln(Number of Groups)"))  + ggtitle("") + 
  scale_colour_manual(values = c(1, 2), name = "Shared Ties") + 
  scale_linetype_manual(values = c(1, 2), name = "Shared Ties") +  
  theme(plot.title = element_text(hjust = 0.5)) +
  theme(legend.position = c(0.7, 0.8))
a
pdf("prednum.pdf", 4, 3)
plot(p_1000, colors="bw") + 
  theme_bw() + ylab("Predicted Pr(External Support)") + xlab(expression("Ln(Number of Groups)"))  + ggtitle("") + 
  scale_colour_manual(values = c(1, 2), name = "Shared Ties") + 
  scale_linetype_manual(values = c(1, 2), name = "Shared Ties") +  
  theme(plot.title = element_text(hjust = 0.5)) +
  theme(legend.position = c(0.7, 0.8))
dev.off()

margins(reg.m1)


write.csv(df2, "externalsupport_robustnessdata.csv")

df2=read.csv("/Users/irismalone/Dropbox/AGDPapers/IS-ISSS Paper/Data/v4/externalsupport_robustnessdata.csv")
df2 = with(df2, data.frame(ccode1,
ccode2,
year,
dyadid,
torgid,
polopp,
anystatesup,
statetypesupmaterial,
statetypesuptrain,
statetypesupfinancial,
statetypesupother,
statetypesupterr,statetypesupendorse,
sharedethnicpower,
leftideational,
sunniideational,
shiaideational,
l1noviolsupl,
l1nostatesupl,
sharedcoideo,
sharedcoideo2,
sharedrelig,
coldwar,post911,
rivalry,
minwbgdppc, 
mindemest,
lndist,
cowmaj2,
aim2,
polwing,
transnationalattackany,
civwar, region.x))
write.csv(df2, "/users/irismalone/Dropbox/AGDPapers/IS-ISSS Paper/Replication Files - II/Data/replication_extsup_data.csv")